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ABSTRACT 

Stars bom from the same molecular cloud should be nearly homogeneous in their element abundances. The 
concept of chemical tagging is to identify members of disrupted clusters by their clustering in element abun¬ 
dance space. Chemical tagging requires large samples of stars with precise abundances for many individual 
elements. With uncertainties of cr^x/vt] ^^d cr[Fe/H] — 0-05 for 10 elements measured for > 10"^ stars, the 
APOGEE DR12 spectra may be the first well-suited data set to put this idea into practice. We find that even 
APOGEE data offer only ^500 independent volume elements in the 10-dimensional abundance space, when 
we focus on the G-enhanced Galactic disk. We develop and apply a new algorithm to search for chemically 
homogeneous sets of stars against a dominant background. By injecting star clusters into the APOGEE data 
set we show that chemically homogeneous clusters with masses ^ 3 x 10^ Mq would be easily detectable and 
yet no such signal is seen in the data. By generalizing this approach, we put a first abundance-based constraint 
on the cluster mass function for the old disk stars in the Milky Way. 

Keywords: Galaxy: abundances — Galaxy: disk — Galaxy: evolution — Galaxy: formation — ISM: abun¬ 
dances — stars: abundances 


1. introduction 

The Milky Way offers a unique opportunity to understand 
how disk galaxies form, in particular when and where, and 
in which types of aggregates, or clusters, they formed their 
stars. As star cluster masses depend on gravitational instabil¬ 
ities (e.g., Escala & Larson 2008), the cluster mass function 
(CME) indirectly probes the dynamical state of the Milky Way 
disk over cosmic time. 

At least during the intensely star-forming, early phases of 
the Milky Way, the majority of stars are believed to form in 
clusters (Kruijssen 2012; Adamo et al. 2015). Most of these 
are rapidly disrupted and dispersed throughout the Galaxy 
(e.g., Odenkirchen et al. 2003; Koposov et al. 2010; Dalessan- 
dro et al. 2015), for a brief while appearing as moving groups 
(e.g., De Silva et al. 2007a, 2013; Bubar & King 2010). Once 
the phase space information as a common birth marker is 
lost, chemical tagging, first proposed by Ereeman & Bland- 
Hawthorn (2002), may still betray the common birth origin 
of stars through their exceptional similarity in element abun¬ 
dances. Stars originating from the same star cluster are be¬ 
lieved to be homogeneous in their chemistry (e.g., De Silva 
et al. 2007b; Koposov et al. 2008; De Silva et al. 2009; Ting 
et al. 2012b; Eeng & Krumholz 2014; Eriel et al. 2014). Since 
the photospheric element abundances, at least for elements 
heavier than sodium, are invariant throughout their lifetime, 
they are permanent tags of the stellar birth origins. 

A broad goal of chemical tagging is to reconstruct the stel¬ 
lar CME, i.e. the relative distribution of (chemically homo¬ 
geneous) stellar cluster masses when they formed. Although 
we can investigate the present-day CME through young star 
clusters and massive surviving clusters (e.g., Bica et al. 2003; 
Lada & Lada 2003; Porras et al. 2003; Brandner et al. 2008; 
Borissova et al. 2011; Bragaglia et al. 2012), the Milky Way’s 
CME in the past is unknown. The key idea in chemical tag- 
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ging is that massive chemically homogeneous clusters show 
up as discernible clumps in the multi-dimensional abundance 
space (e.g., Bland-Hawthorn et al. 2010a; Ting et al. 2015). 

Besides understanding the CME, chemical tagging is 
deemed an essential tool to understand the effect of radial mi¬ 
gration in Galactic evolution (e.g., Schonrich & Binney 2009; 
Gould & Rix 2015). Sellwood & Binney (2002) first proposed 
that stars could migrate significantly from their radial position 
when resonate with spiral/bar structures (see observational ev¬ 
idence from Haywood 2008; Loebman et al. 2011; Kordopatis 
et al. 2015). Although simulations concur to the analytic cal¬ 
culations (Roskar et al. 2008, 2012; Bland-Hawthorn et al. 
2010b; Minchev & Eamaey 2010; Di Matteo et al. 2013; Halle 
et al. 2015), direct observational evidence of radial migration 
remains controversial. The ability to recover dispersed star 
clusters would be fundamental in quantitatively constraining 
radial migration models. 

In recent years, the idea of chemical tagging has garnered 
more attention. Large spectroscopic surveys, including RAVE 
(Steinmetz et al. 2006), APOGEE (Zasowski et al. 2013), 
GALAH (De Silva et al. 2015) and Gaia-ESO (Randich et al. 
2013), are being carried out. Results from these surveys 
have demonstrated the power of using elemental abundance 
patterns to identify distinct stellar populations in the Milky 
Way. Martig et al. (2015) found young G-enhanced stars that 
are difficult to explain within current models of the evolu¬ 
tion of the Milky Way. Masseron & Gilmore (2015) and 
Hayden et al. (2015) showed that the C/N ratio and the g- 
elements are good indicators to separate the thin and thick 
components of the Galactic disk. Schiavon et al. (submit¬ 
ted) found bulge stars that show abundance patterns similar 
to globular clusters. These studies focus on finding popula¬ 
tions of stars via their abundance patterns, which is a “weak” 
form of chemical tagging. The goal of the “strong” form of 
chemical tagging is to identify stars that were born from the 
same molecular cloud. When we refer to chemical tagging in 
this paper, we only refer to this “strong” form. In this con¬ 
text, due to these exciting new opportunities, many prepara- 
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tory works started to explore the capabilities of these surveys 
and the idea of chemical tagging. For example, Tabernero 
et al. (2012, 2014); Mitschang et al. (2013, 2014); Blanco- 
Cuaresma et al. (2015); Macfarlane et al. (2015) proposed var¬ 
ious schemes and performed numerical experiments on sepa¬ 
rating open clusters/moving groups in abundance space. 

Ting et al. (2015) explored the feasibility of chemical tag¬ 
ging by exploring a grid of Galactic evolution parameters. 
They found that identifying individual clusters through chem¬ 
ical tagging is generally challenging with the on-going sur¬ 
veys. Clumps that show overdensities in abundance space are 
usually made up of many smaller clusters, i.e., the background 
contaminants are non-negligible in clump search. One might 
be able to associate a detected clump as a single disrupted 
cluster only when the background density in abundance space 
is low. They showed that a low background can be achieved 
by either studying a subpopulation that occupies a large vol¬ 
ume in abundance space, small abundance uncertainties or 
a large number of independent elements. Despite all these 
complications to reconstruct individual clusters, they argued 
that we can still statistically reconstruct the CMF through the 
dumpiness in abundance space. The main goal of this paper 
is to apply that idea to the APOGEE data (Holtzman et al. 
2015). 

In §2, we characterize the APOGEE sample that we explore 
in this study. In §3, we introduce our clump search method. 
The key is to define a robust search sphere that has the high¬ 
est signal-to-background ratio possible. In §4, we apply this 
method to the APOGEE DR 12 data. We discuss how this re¬ 
sult can be applied to obtain a tentative constraint on CME. 
We conclude in §5. By comparing the signal-to-background 
contrast observed in the data to the simulated contrast from 
injected clusters, we will argue that no chemically homoge¬ 
neous clusters more massive than 3 x 10^ Mq have formed in 
the Qf-enhanced disk. 

2. APOGEE SAMPLE PROPERTIES 

We adopt the APOGEE DR12 publicly available sample 
(Holtzman et al. 2015). Similar to Hayden et al. (2015), we 
consider stars with all elements measured and with reliable 
abundances, i.e., 4,000 < Teff < 5,500, 1 < logg < 3.8 and 
signal to noise ratio > 80. Stars from APOGEE that satisfy 
these criteria are plotted in Eig. 1. In this study, we only fo¬ 
cus on the (T-enhanced disk as defined via the cut shown in 
Eig. 1. We focus on this subsample as their chemical/spatial 
modeling is likely to be more straightforward (see §4.2 and 
also see Bland-Hawthorn et al. 2014; Ting et al. 2015). Eur- 
thermore, the Milky Way was likely kinematically hotter and 
more turbulent in the first few billion years (e.g., Boumaud 
et al. 2009; Krumholz et al. 2012; Bird et al. 2013). High 
redshift (z ^ 2) extragalactic studies have revealed the exis¬ 
tence of massive star-forming clumps in star-forming galax¬ 
ies (e.g., Livermore et al. 2012; Genzel et al. 2013). As a 
result, star clusters within the ^-enhanced disk could be more 
massive and if so would provide a strong signal in abundance 
space. The ^-enhanced disk also occupies a larger volume 
in abundance space, i.e., lower background density of stars, 
which guarantees clumps a better contrast to the background 
in abundance space. 

We checked that our main result presented in this study, 
namely there is no cluster more massive than 3 x 10^ Mq 
formed in the Milky Way, still holds at least to the Icr level 
(see §4.3.2), if we choose a selection cut within the shaded 
region as shown in Eig. 1. The result only changes more dra- 



Figure 1. Black symbols show 13,000 APOGEE stars selected for this 
study. We only consider a-enhanced stars, as the volume they occupy in 
10-dimensional abundance space is larger than the low-a sequence, so the 
background will be lower, and therefore detection is more likely. The red 
symbols show the ~ 7% of outliers in the 10-dimensional abundance space 
that we do not include in the sample. After we discard these outliers, the 
10-dimensional empirical distribution in abundance space is better modeled 
by an ellipsoid and is easier to deconvolve (see §3.2). The solid black line 
shows the fiducial selection cut. We also examine that the results in this paper 
are not sensitive to our data selection. If we choose a selection cut within the 
shaded region, the results in this paper still hold. The results only change 
more dramatically if we choose a much steeper cut such as the blue dashed 
line (see text for details). 

matically if we choose a much steeper cut such as the blue 
dashed line. In this case, we discard too many low density 
regions where most chemical tagging signals reside (see §3.1 
and §4.3), and we can only rule out clusters > 10^ Mq. Using 
the fiducial selection cut, in total, the ct-enhanced sample has 
14,002 stars, as shown in the red and black symbols in Eig. 1. 
Eor reasons and selection criteria that will become clear in 
§3.2, we further discard the 7% of most outlying stars, shown 
as red symbols in Eig. 1, and end up with a final sample of 
13,000 stars. 

Without further kinematic modeling, it is hard to disentan¬ 
gle the halo stars from the disk stars. But the elimination of 
outliers as shown in Eig. 1 culls most of the metal-poor stars 
with [Ee/H] < -1 and therefore, the bulk population in this 
study should not be contaminated much by the halo stars. We 
also performed the same analysis only considering stars with 
[Ee/H] > -1 and found that the results remain qualitatively 
the same. Eor the potential bulge contamination, we find that 
among the 13,000 stars, only 3% of them satisfy the bulge 
stars criteria with I <22, |Z?| < 15 and Rgc < 3kpc, where 
Rgc is the isochrone Galactic radius derived in Hayden et al. 
(2015). Therefore, we will assume throughout this study that 
the sample only consists of disk stars. 

To perform the chemical tagging experiment, we want to 
consider as many elements as possible in order to maximize 
the volume in abundance space. In this case the background 
becomes more diluted, and the signals will therefore have a 
better chance of standing out from the background. In to¬ 
tal, APOGEE measures 15 elements. However, as discussed 
in Holtzman et al. (2015), Na, Ti, V abundances might not 
be reliable yet in the current release. Eurthermore, C and N 
are expected and seen to evolve through stellar evolution due 
to the post main sequence dredge up (e.g. Karakas & Lat- 
tanzio 2014; Masseron & Gilmore 2015; Ventura et al. 2015). 
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Therefore, it is complicated to relate C and N to their prim¬ 
itive abundances when the stars formed. Discarding these 5 
elements, we consider 10 elements in this study, namely Al, 
K, O, S, Mg, Si, Ca, Mn, Ni, Fe. Next, we want to define 
an abundance space out of these 10 elements. As discussed 
in Ting et al. (2012a), chemical tagging is better performed 
in [X/Fe] instead of [X/H]. [X/H] strongly correlate with each 
other. It is harder to observe the subtle variants among clus¬ 
ters in an abundance space spanned by [X/H]. Therefore, in 
this study, we consider an abundance space spanned by [Fe/H] 
and 9 [X/Fe] from elements beside iron. We denote a vector 
in this 10-dimensional space to be X in this study. 

As we will discuss in more detail in §3, to find a chemi¬ 
cally homogeneous cluster in abundance space we first need 
to understand the typical volume that such a cluster occupies, 
after accounting for the measurement uncertainties that will 
dominate over the intrinsic scatter. Therefore, to estimate the 
typical volume, we will evaluate differential uncertainties, or 
measurement precision (not accuracy), ax and their correla¬ 
tions, i.e., the empirical point spread function of a chemically 
homogeneous cluster in abundance space. We estimate this 
co-variance matrix from known clusters in the APOGEE data 
and refer to the resulting matrix as the “cluster kernel” below. 

We consider the DR 10 clusters classification (Meszaros 
et al. 2013) since the DR 12 classification is yet to be re¬ 
leased. We cross-match the DRIO cluster member identi¬ 
ties with DR 12 and adopt the element abundances from the 
DR 12 release. We only consider clusters with more than 10 
cluster members. Three open clusters (NGC6819, NGC2158, 
M67) satisfy this criterion. Noting that all these clusters are 
metal-rich with [Ee/H] ^-0.1 and the possibility that abun¬ 
dance determination could be worse at lower metallicity, we 
also adopt one of the more metal-rich globular clusters. Ml07 
([Ee/H] -1). We fit Ml07 members with two Gaussians 

distributions in the 10-dimensional abundance space to elim¬ 
inate any possible secondary population from this globular 
cluster. Including Ml07, we have a total of 77 cluster mem¬ 
bers. The primary population of Ml07 shows measurement 
uncertainties consistent with other open clusters. Restrict¬ 
ing ourselves to the three metal-rich open clusters results in 
a slightly smaller ax- Therefore, including Ml07 only makes 
our results more conservative (see §4.3.1). We subtract the 
element abundances of each cluster by their means to center 
clusters at the zero origin. The co-variance matrix of these 
77 stars is estimated. This matrix defines an ellipsoid that a 
typical chemically homogeneous cluster occupies. Eor each 
element, we find that ax ^ 0.05-0.06 dex; this multivariate 
Gaussian sets the effective volume that homogeneous clus¬ 
ters occupy in the observed abundance space. This estimate 
is consistent with Holtzman et al. (2015) (see table 6 in the 
paper, but note that they show measurement uncertainties in 
[X/H], but here we evaluate uncertainties in [X/Ee] [Ee/H]). 
Due to the small sample size of cluster members, we boot¬ 
strap the cluster sample and find that the uncertainty on this 
ax estimate is about 20%. A larger sample of cluster calibra¬ 
tors would be very helpful as a precise estimate of the cluster 
kernel is a key ingredient in any chemical tagging analysis. 

3. METHOD 

In this section, we will describe the challenges in abundance 
clump searches and our clump search method. Although 
various schemes have been proposed to separate open clus¬ 
ters/moving groups (e.g., Sharma & Johnston 2009; Blanco- 
Cuaresma et al. 2015; Macfarlane et al. 2015) in abundance 



Figure 2. Distributions of all 10 elements in this study. By comparing these 
distributions to the effective search sphere diameter (see §3.1 for details), we 
find the abundance distributions in each dimension to have standard devia¬ 
tion widths that are typically ~ 1 - 3 times of the effective search sphere ra¬ 
dius. The relatively small volume in abundance space compared to the search 
sphere volume highlights the main challenge in chemical tagging. 

space, a question often not discussed is the estimation and in¬ 
clusion of background contaminants. Simulations from Ting 
et al. (2015) showed that the background contaminants can 
be a critical limiting factor in chemical tagging experiments. 
After extensive experimentation we found that most proposed 
techniques, such as K-means, Gaussian mixture models, and 
minimal spanning tree, are only effective in separating clumps 
in the limit of a small background or when the background can 
be easily estimated. 

Due to this limitation, we have developed a simple new 
method^ geared toward regimes where the background is 
dominant and has a complex topology in abundance space 
(read also Everitt et al. 2010). The key to our method consists 
of two parts that we will explain in §3.1 and §3.3. Eirst, we 
need to estimate the local density. As we will discuss in more 
detail in §3.3, we define the local density to be the number of 
stars within a search sphere. The search sphere that we use to 
distinguish signals from the background should be sufficiently 
large. It should include a large fraction of a chemically ho¬ 
mogeneous cluster but avoid being too wide and should not 
include too many background contaminants. Secondly, the 
abundance space distribution is not uniformly distributed. To 
estimate the detection significance, we have to estimate the 
expected background at each location. 

3.1. Abundance space search sphere 

In §2, we derived the empirical multivariate Gaussian dis¬ 
tribution that a typical cluster occupies in abundance space. 
This distribution defines an ellipsoidal distribution in the 10- 
dimensional abundance space. Since an optimal search sphere 
should include a high ratio of cluster objects compared to 
the background contaminants, the optimal search volume, or 
search sphere, should follow the same uncertainty ellipsoid. 
However, it is not convenient to operate with a tilted ellipsoid 
because a simple Euclidean distance from the center alone is 
not sufficient to determine whether a star is included in the el¬ 
lipsoid. Therefore, we linearly transform the abundance space 

^ Our method can be regarded as a variation of density-based nonparamet- 
ric clustering techniques. 
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such that the cluster kernel becomes a unit Gaussian distribu¬ 
tion. We emphasize that the transformation is only to make 
calculations more straightforward, leaving all astrophysical 
implications invariant. 

We can now determine an appropriate radius for the search 
sphere in these coordinates. A unit radius is not a good choice 
for the search sphere even though clusters follow a unit Gaus¬ 
sian distribution in the transformed coordinates: a box with 2 
dex in width in each dimension only captures (68%)^^ = 2% 
of the clump. Since a unit n-sphere is strictly included in this 
box, one can show that a unit n-sphere encapsulates an even 
smaller fraction, 0.02%, of the data. This curse of dimen¬ 
sionality"^ implies that to capture, for example, 68% or la of 
the cluster members, we require a search sphere with a ra¬ 
dius larger than 1 dex in the transformed abundance space. In 
fact, mathematically, the inclusion fraction of a unit Gaussian 
within an n-sphere follows the -distribution. The analytic 
formula of a -distribution shows that an n-sphere of 3.4dex 
in radius is needed. 

One important parameter that will determine the difficulty 
of chemical tagging detections is the number of separate 
“chemical cells” in the abundance space (Freeman & Bland- 
Hawthorn 2002; Ting et al. 2015). The number of chemi¬ 
cal cells is the ratio between the abundance space volume 
spanned by the data and the typical volume of a search sphere. 
We will evaluate this number later. But for now, one way to 
visualize the number of chemical cells is to compare the ef¬ 
fective diameter of the search sphere to the spread in each 
element abundance. We illustrate this comparison in Fig. 2. 
We find that an n-sphere with a radius of 3.4dex in the trans¬ 
formed abundance space corresponds to an ellipsoid with an 
effective radius of 0.05 dex in the original abundance space. 
The effective radius is defined such that an n-sphere with this 
radius contains the same volume as the ellipsoid. Due to the 
large dimensionality, we find that if we do not take into ac¬ 
count the covariances of the cluster kernel, i.e., if we were to 
use an n-sphere in the original abundance space, instead of a 
tilted ellipsoid defined from the cluster members, we estimate 
that the search volume would be 10,000 times larger and the 
search sphere would include too many background contami¬ 
nants. 

Even with this optimized search ellipsoid in the original 
abundance space, as shown in Fig. 2, the distribution of each 
element is typically only 1-3 times the effective search sphere 
diameter. As a result, it is not possible to search for clusters 
in the core region of the chemical distribution. In this region, 
the search sphere would include too many background con¬ 
taminants. The chemical tagging signals are most likely to 
come from the peripheral regions of the chemical distribution 
where the background contaminants are not dominant (also 
see Bland-Hawthorn et al. 2010a; Karlsson et al. 2012; Bland- 
Hawthorn et al. 2014). Nonetheless, in a 10-dimensional 
space, the “surface-to-core” ratio is very large. Therefore, 
there is a reasonable chance of finding clumps in the periph¬ 
eral regions. As we will show in §4, all chemical tagging 
signals indeed come from the peripheral regions. 


^ Techniques to compactify of dimensions, such as PCA (e.g., Ting et al. 
2012a), do not mitigate this problem because the density of background con¬ 
taminants also increases accordingly in the compactified space. 


3.2. The intrinsic abundance distribution of the a-enhanced 

sample 

We now proceed to deconvolve the observed abundance 
distribution. The intrinsic abundance distribution is required 
in order to inject mock clusters into the observed APOGEE 
dataset. 

In the previous section, we transformed abundance space 
such that clusters follow a 10-dimensional unit Gaussian dis¬ 
tribution. However, in such coordinates, the overall chem¬ 
ical distribution of the ^-enhanced disk will still presum¬ 
ably show co-variances among the different coordinate di¬ 
rections because the transformation is only to normalize the 
cluster point spread function and makes no assumption on 
the global distribution. Deconvolving such a co-variant dis¬ 
tribution directly in 10-dimensions is computationally pro¬ 
hibitive. Therefore, we further rotate the transformed abun¬ 
dance coordinate system to eliminate the co-variances such 
that the joint 10-dimensional abundance distribution of the 
(T-enhanced disk can be approximated by a product of 10 
marginal distributions. Since the cluster kernel, refiecting the 
measurement uncertainties, is already isotropic in the trans¬ 
formed abundance space, it remains unaffected by further ro¬ 
tation. Upon this rotation, the deconvolution task simplifies 
to ten independent, one-dimensional deconvolutions on the 
marginal distributions. 

However, this approach only works if the actual ensem¬ 
ble abundance distribution can be well approximated by its 
marginal distributions after rotation, i.e. if it does not show 
significant curvature in abundance space. If we consider all 
14,002 (T-enhanced stars in the sample, we find that this is 
not a good approximation: the 10-dimensional abundance 
space has a less regular topology, dominated by a small frac¬ 
tion of outliers in abundance space, as shown in red sym¬ 
bols in Fig. 1. To look for these outliers, we perform a 10- 
dimensional kernel density estimation, using a unit Gaussian 
distribution as the smoothing kernel. We rank the data points 
according to their local density in the kernel density estima¬ 
tion and discard the most outlying 1,002 stars. We check that 
upon discarding these outliers (7%), injecting clusters accord¬ 
ing to the joint distribution gives similar statistics to injecting 
clusters according to product of marginal distributions. We 
emphasize that discarding outliers shrinks the peripheral re¬ 
gions and makes detection more unlikely. The main purpose 
of this study is to put an upper limit on the maximum cluster 
mass (see §4.3.2), discarding outliers only makes our estimate 
more conservative. We also check that these outliers are not 
particularly clumped in abundance space and hence are un¬ 
likely to be chemical tagging detections. 

After breaking down the joint distribution to its marginal 
distributions, ^convolved,/, we model each marginal distribution 
with the sum of two Gaussians (see also McLachlan & Peel 
2000; Everitt et al. 2010), 

^convolved, 

~ {\-f)M{Xi\ii\,a\)+fN{Xi\n2,o-l). ( 1 ) 

where p and a are the means and standard deviations of the 
Gaussian distributions, / shows the relative contribution from 
each distribution. We require a two-component Gaussian be¬ 
cause the marginal distributions often show a core region and 
a broad wing region. Fitting a single Gaussian will underesti¬ 
mate the total area of low-density wings. As discussed in §3.1, 
the wing regions are the most valuable parts of a chemical dis¬ 
tribution as they have the highest chance to detect clumps. In 
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Figure 3. Deconvolved (intrinsic) model of the 10 elements in this study, projected onto [A/Fe]-[Fe/H] planes. The black symbols show the APOGEE sample 
adopted, illustrating the empirical distribution. The red contours show the 50 and 90 percentiles of the projected intrinsic distributions. When injecting mock 
clusters, we draw their center locations from the intrinsic distribution model instead of the empirical APOGEE distribution. Note that although an ellipsoid model 
provides an acceptable fit to the data, in some cases the abundance distribution is influenced by non-ellipsoidal structures. 


the model, we allow centers of the two Gaussian distributions, 
/ii and /i 2 , to be different to account for an asymmetric dis¬ 
tribution. We found that this model provides an excellent fit 
to each ID marginalized distribution (although the joint dis¬ 
tribution fits are slightly less satisfactory as we will explain 
below) and the deconvolution can be done analytically. We 
model the intrinsic chemical distribution /^intrinsic to be 

10 

^intrinsic(^) “ ^intrinsic ,/(-^0 5 ( 2 ^) 

dim=l 


where 

^intrinsic,/ (-T/) = ^convolved,/ (-T/1/^1, /^2 1 ,(J2~1 (3) 

Fig. 3 shows the intrinsic chemical distribution derived ac¬ 
cording to the method above. We will use this model to 
draw mock data of hypothetical clusters, whose abundance 
probability is drawn from the ensemble distribution. We cau¬ 
tion that the 10-dimensional ellipsoidal model does not give 
a perfect fit to the data, despite the fact that we have elimi¬ 
nated 7% of the outliers. For instance the [Ca/Fe] vs. [Fe/H] 
distribution, as illustrated in Fig. 3, shows a more complex 
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morphology than an ellipsoidal model. In the ideal case, we 
would draw hypothetical clusters from a deconvolved distri¬ 
bution that displays similar intricate morphology. However, 
deconvolving such intricate distribution is computationally in¬ 
tractable in 10-dimensional. Nonetheless, we checked that in¬ 
jecting clusters according to the convolved ellipsoidal model 
shows similar statistics to the empirical convolved distribu¬ 
tion. Therefore, in this study, we make the assumption that in¬ 
jecting cluster according to the deconvolved ellipsoidal model 
gives similar statistics if we were to inject according to the 
real deconvolved distribution. 

3.3. Detection significance 

So far, we have defined an operative search sphere with ra¬ 
dius r = 3.4dex to look for overdensities in the transformed 
abundance space. We know that the overall abundance distri¬ 
bution, i.e. the background, is not uniform (e.g., Edvardsson 
et al. 1993; Barklem et al. 2005; Reddy et al. 2006; Bensby 
et al. 2014). The absolute number of stars within the search 
sphere is therefore not particularly informative. We want to 
find regions where the local density within a search volume 
is significantly higher than its vicinity regions. Therefore, to 
define the detection significance, we need two ingredients: a 
local density estimation at each location and the correspond¬ 
ing local background estimation at this location. 

Fig. 4 shows a schematic illustration of our clump search 
method. For each star 5', we define the local density of this 
star, ns, to be the total number of stars located within r dis¬ 
tance from this star. Throughout this study, we only consider 
stars with > 10 to avoid the large Poisson fiuctuation at 
small ns. We estimate the vicinity background, (ns) to be the 
average of ns' where s' are all stars located within a distance 
of r-2r from star 5'. We define the detection significance to be 
c^detect(‘5’) = (ns-(ns))/\J(ns), which measures the deviation of 
the local density from the vicinity background, in units of the 
Poisson uncertainty of the background. 

Although this is a sensible definition of detection signifi¬ 
cance, there is a complication. If we have a uniform back¬ 
ground, provided there is no signal, (idetect should center 
around zero. Unfortunately, this is not the case for an un¬ 
even background, especially for high dimensions. At a fixed 
point in an uneven background, there are always more vicin¬ 
ity regions that have lower densities (toward the valley) than 
regions that have higher densities (toward the core). As a re¬ 
sult, we have ns > in general. This disproportion gets 
more severe toward the core as the background gradient gets 
steeper. This disproportion causes crdense(^ 5 ) to be an increas¬ 
ing function of ns. To overcome this shortcoming and to have 
c^detect centered around zero, we calibrate (idetect by the median 
of crdetect(^ 5 ) at each ns. We denote the calibrated detection 
significance to be ^detect and use it to be our operative mea¬ 
sure of detection significance in the following and apply this 
method to the APOGEE data. 

4. RESULTS 

We now explore what we can learn about the number of 
chemical cells in the APOGEE survey, and about the presence 
of any clumps in abundance space, that may reflect chemically 
tagged remnants of disrupted clusters. 

4.1. The Number of chemical cells in the APOGEE 
observations of the a-enhanced Galactic disk 

The number of chemical cells is best estimated in the trans¬ 
formed (and rotated) coordinates where the global chemical 


Star / Location 



Stars 
in vicinity 


Figure 4. Schematic illustration of our clump search method. At each star 
5 (white symbol), we evaluate the number of stars (orange symbols) within 
the search sphere (red shaded region) that we denote as Ug. The vicinity back¬ 
ground, (n^/} is calculated by averaging other , as shown in the blue shaded 
regions, where s' are all stars (yellow symbols) that are in the vicinity region 
(green shaded region). The vicinity region is defined to be the region out¬ 
side the search sphere but inside two times the search sphere. The detection 
significance is then defined as crdetect(‘^) = ))/ 

distribution has no co-variances between different coordinates 
and a chemically homogeneous cluster can be represented 
by a unit Gaussian distribution. To calculate the number of 
chemical cells, let us estimate the global distribution to be a 
multivariate Gaussian with no co-variances and with standard 
deviations cri,cr 2 ,... ,crio in the 10 transformed coordinates. 
The number of chemical cells, by definition, is the volume 
ratio between the global distribution over the cluster kernel. 
Since the cluster kernel has a unit width in all directions, the 
number of chemical cells in APOGEE can be estimated to be 
(cri-cr2---CTio)/(ldex)“~500. 

This estimate agrees with the prediction in Ting et al. 
(2012a) using principle components analysis. The lack of 
chemical cells despite having 10-dimensions is due the strong 
correlations among abundances, especially for the a-capture 
elements and Fe peak elements. The small number of chem¬ 
ical cells emphasizes the challenges in performing chem¬ 
ical tagging with strongly correlated element abundances. 
Nonetheless, as we have discussed in §3.2, the APOGEE 
abundance space has broad wings that are not captured in a 
single multivariate Gaussian. Therefore, in the analysis pre¬ 
sented below, the signals drawn from a composite multivari¬ 
ate Gaussian, as described in equation (1), are stronger than a 
simple multivariate Gaussian with 500 chemical cells. 

4.2. Relation between A^mject nnd Mduster 

To understand whether a cluster will be detected or not, we 
first need to investigate the number of stars that we would 
sample in APOGEE from a cluster, given its zero age clus¬ 
ter mass Mduster- We denote the number of stars sampled to 
be Mnject and will use the one-to-one relation between Ainject 
and Mduster in the following discussion. But this conversion is 
based on some critical assumptions that we will now explain. 
The relation between Ainject and Mduster is one-to-one up to a 
Poisson uncertainty - more massive clusters have more stars 
to begin with, and therefore will have more stars sampled in 
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Figure 5. Calibrated detection significance as a function of local density in abundance space. The left panel shows the sample observed values. The dashed lines 
show the observation detection baseline. The overlaid orange symbols in the middle and right panels show results from the injected 10^ Mq and 10^ Mq clusters, 
respectively. Since we calculate the local density centered at each star, there are numerous points for each cluster. We highlight the maximum deviation of each 
cluster with a bold red symbol. We compile results from 1000 and 100 trials for the middle and right panels. To demonstrate a typical detection significance 
distribution within a cluster, we link with solid gray lines all stars from the cluster containing the star with the highest detection significance. The middle and 
right panels show that if a 10^ -10^ Mq cluster formed in the past, some cluster members might show more deviations than the detection baseline. 


the survey. 

In the limit where there is no radial migration, the relation 
between A^mject and Mduster is simple and can be derived analyt¬ 
ically. Assuming stars are azimuthally mixed in the annulus, 
the number of stars sampled from a cluster, A^mject. can be ap¬ 
proximated (see Bland-Hawthorn et al. 2014; De Silva et al. 
2015; Ting et al. 2015) to be 


A^i 


inject ■ 


Afcluster 
Afannulus 


A^apogee, 


( 4 ) 


where Mannuius is the total stellar mass (including stellar mass 
loss) in the annulus and A^^pqgee = 13,000 is the APOGEE 
sample size in this study. This formula can be easily under¬ 
stood as the following. Assuming stars in the sample have an 
average stellar mass (M) IMq, the ratio A^APOGEE/^^annuius 
gives the stellar mass fraction within the annulus that we 
would sample in the survey. We denote this ratio to be the 
sampling rate. When multiplying the sampling rate by a clus¬ 
ter mass, the product gives the stellar mass, and thus the num¬ 
ber of stars with (M) = 1 Mq, that we would sample from this 
cluster. 

However, radial migration modifies Mannuius in n complex 
way (see details in Ting et al. 2015). Stars born outside the 
annulus could migrate into the annulus, and stars from the an¬ 
nulus could now appear to be outside the annulus. Due to 
this complication, to estimate Mannuius, the APOGEE’S selec¬ 
tion function, as well as a robust Galactic chemical evolu¬ 
tion model (e.g., Kobayashi et al. 2006, 2011; Minchev et al. 
2013), is needed. This is clearly beyond the scope of this pa¬ 
per. 

To simplify the problem and to only derive a conservative 
limit on the CME in §4.3.2, we assume that the a-enhanced 
disk is completely radially mixed. In other words, the current 
spatial location of a star is completely random and is inde¬ 
pendent of their birth radii. In this limit, a star in the sam¬ 
ple can be any star from the a-enhanced disk. Therefore, we 
have Mannuius = Mtotai, where Mtotai is the total stellar mass of 
the (T-enhanced disk. Although complete mixing is a crude 
assumption, it is likely to be reasonable for the ^-enhanced 
disk. Eor example, Hayden et al. (2015) showed that there is 
a universal a-trend irrespective of the Galactocentric radii. A 
natural explanation of this result is the stars in the a-enhanced 
disk are well-mixed. Moreover, the APOGEE sample covers a 


wide range of Galactocentric radius, with Rgc — 3kpc- 15kpc 
(Bovy et al. 2014; Nidever et al. 2014; Hayden et al. 2015). It 
should have sampled the ^-enhanced disk from a large frac¬ 
tion of the Milky Way. 

We emphasize that the complete mixing assumption gives 
a conservative limit on the CME. In the case where the mix¬ 
ing is not complete, we would have sampled more stars from 
the same cluster, and hence, it would be easier to exceed the 
APOGEE baseline. On top of the complete mixing assump¬ 
tion, we also assume that the CME is independent of element 
abundances, and hence the equation (4) applies universally to 
the whole abundance space. 

With these assumptions, we only need to properly estimate 
Mtotai and apply equation (4) to obtain a one-to-one relation 
between Abject and Mduster- We assume that the ct-enhanced 
disk has an exponential scale length of 3kpc (e.g., Bovy et al. 
2012) and consists of 10% stellar mass observed in the solar 
neighborhood (e.g., Cheng et al. 2012). We adopt the stel¬ 
lar density in the solar neighborhood to be 38 MqPC“^ (e.g., 
Elynn et al. 2006; Bovy & Rix 2013; Zhang et al. 2013), and 
the solar Galactocentric radius Rq = 8kpc (e.g., Ghez et al. 
2008; Gillessen et al. 2009; Reid et al. 2014). These as¬ 
sumptions yield a present-day ^-enhanced disk stellar mass 
of ^ 3 X 10^ Mq. Since the ^-enhanced disk is old, massive 
stars have long since evolved and died. To account for this, we 
consider a total stellar mass loss of 40% (Conroy et al. 2009, 
assuming a Kroupa IME). Putting all these together, we have 
Mtotai — 6 X 10^ Mq. We also derive the sampling rate of the 
current APOGEE sample to be Aapogee /Mtotai = 5 ^\q 5 • On 

average, we would collect one star from a 5 x 10^ Mq cluster. 
We will defer the discussion on what this low sampling rate 
implies in §4.3.2. 

4.3. Chemical tagging in APOGEE 

We apply the clump search method described in §3 to the 
APOGEE sample. The left panel of Eig. 5 shows ^detect as 
a function of ns of all 13,000 stars. At face value, it is tan¬ 
talizing to observe deviations > 5cr. But we emphasize that 
the detection significance depends on the various assumptions 
made, such as the choice of search sphere radius, the mini¬ 
mum number of neighbors requirement, the detection signifi¬ 
cance calibration and the definition of vicinity region at each 
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data point. Without further information such as stellar ages, it 
is difficult to confirm the origin of these clumps. Furthermore, 
as shown in Ting et al. (2015), most clumps are comprised of 
many clusters sharing similar element abundances. 

Due to these uncertainties, instead of interpreting these 
clumps as detections, we proceed by assuming the APOGEE 
dataset (left panel of Fig. 5) to be the detection baseline. We 
inject mock clusters of different sizes into the data and esti¬ 
mate the detection significance of these injected objects. By 
forward modeling, we rule out cases that are not consistent 
with the observation baseline. The middle and right panels 
show ^detect of the injected objects. In the right panel, we 
combine results of 100 trials, where in each case we inject a 
IO^Mq (Mnject — 250 stars) clump into the data. In the middle 
panel, we show the results of 1000 trials with IO^Mq clus¬ 
ters (Mnject — 25 stars) injected. These two panels show that if 
10^ - IO^Mq clusters have formed in the past, there is a rea¬ 
sonable chance that we would have detected larger deviations 
than the value observed. A cluster with IO^Mq lies above the 
detection boundary about ^ 7% of the time, and a cluster with 
IO^Mq is detected about ^ 30% of the time. 

Not all high mass clusters will exceed the detection base¬ 
line. As shown in the middle and right panels, most clusters, 
especially at high ns, blend into the background. To make ro¬ 
bust statements, we now proceed to quantify the probability 
of a cluster exceeding the observation baseline. 

4.3.1. Detection probability of individual clusters 

In this section, we will quantify the probability of an in¬ 
jected clump exceeding the observation baseline. There are 
two key parameters that determine this probability: (a) the 
number of stars injected as a clump, A^mject, and (b) the clus¬ 
ter location in abundance space. As for the latter, the clump 
centers are drawn from the intrinsic distribution model de¬ 
scribed in §3.2. As for the former, we consider a grid of 
Mnject, ranging from 3 -1000 stars with a step-size of 0.2 
in log scale. When injecting a mock cluster with A^mject 
stars, we allow a Poisson fluctuation of Y^Mnject- The A^mject 
range in this study roughly corresponds to cluster masses of 
IO^Mq -5 X IO^Mq. We run 10"^ trials for each A^mject and 
find that the Monte Carlo uncertainty is negligible with this 
many trials. 

We model the cluster location by ranking all 10"^ trials by 
their Pintrinsic(Xcenter) value. Xcenter IS the clump Center location 
in abundance space. We put the ranking into a linear scale, 
which we will denote as p(Xcenter) ^ [0,100], where 

. _ # (^trials < ^intrinsic(Xcenter)) 

P(^center) = • t x • (^) 

#(trials) 

If the cluster is located near the background dominated core, 
it has a higher value in ^intrinsic because the background den¬ 
sity is very large, and we assign a high p(Xcenter)- Whereas, if 
the cluster is located at the peripheral regions, it has a lower 
^intrinsic value since the background density is low, and we as¬ 
sign a low p(Xcenter). 

For each trial, we inject a clump and estimate the local den¬ 
sity and vicinity background for all objects from the injected 
clump the same way in §3.3. We define a clump to have ex¬ 
ceeded the observational baseline if the maximum detection 
significance of this clump exceeds the baseline as demarcated 
by the dashed lines in Fig. 5. We take the maximum signif¬ 
icance of the whole clump because not all injected objects 


Zero Age Cluster Mass, 



1.0 1.5 2.0 2.5 3.0 

Number of Injected Stars, Log Nj^ject 



Number of Injected Stars, Nj^jg^ 

Figure 6. Probability of an injected cluster showing more significant detec¬ 
tion than the APOGEE data. This probability depends on two parameters: 
the number of injected stars per cluster and the cluster location in abundance 
space. The top panel shows the probability as a function of these two param¬ 
eters. As the number of injected stars increases or the cluster is injected in a 
more peripheral region of abundance space, the chance to exceed the observa¬ 
tion baseline increases. The bottom panel shows the probability marginalized 
over the cluster location, i.e., the probability of detecting a cluster of a cer¬ 
tain cluster mass if the cluster location is randomly drawn from the intrinsic 
abundance distribution. The solid line shows result assuming the best cluster 
kernel estimation, (Jx, as also applied to the top panel. The dashed and dotted 
lines show the marginal probability assuming ±20% statistical uncertainties 
of the ax estimate due to the small sample of cluster members (see §2 for 
details). 


will have high deviations, as shown in the solid gray lines in 
Fig. 5. Objects located near the surface of a clump will blend 
into the background. Only the objects near the clump center 
will have high deviations because the search sphere includes 
a large fraction of the clump. 

The top panel of Fig. 6 shows the probability of exceeding 
the APOGEE baseline. Among all trials that have a similar 


Detection Fraction, 
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cluster location and a similar number of stars, we evaluate the 
fraction of them exceeding the baseline. The x-axis shows 
Mnject and the y-axis shows the cluster location quantified by 
percenter)• This panel illustrates that as the number of injected 
objects increases or the cluster location is increasingly toward 
the peripheral regions, the chance of exceeding the baseline 
improves, consistent with our intuition from Fig. 5. 

The solid line in the bottom panel of Fig. 6 shows the proba¬ 
bility marginalized over the cluster location, i.e., the probabil¬ 
ity of a cluster exceeding the APOGEE baseline as a function 
of its cluster mass if the cluster location is randomly drawn 
from the intrinsic abundance distribution. The marginalized 
probability shows that clusters less massive than 10^ Mq have 
negligible chances of exceeding the baseline, but clusters 
more massive than ^ 10^ Mq begin to show tension with the 
deviations observed in APOGEE. The bottom panel also il¬ 
lustrates that even for a cluster as massive as ^ 5 x 10^ Mq 
( inject — 1000), only about half of the time will a cluster ex¬ 
ceed the baseline. The lack of significant detection from the 
other half is not unexpected. As illustrated in the right panel 
of Eig. 5, if a cluster is located in the core region (i.e., high 
Us), the background becomes dominant. In this regime, most 
objects within the search sphere come from background con¬ 
taminants. Therefore, in the core region, the signal tends to 
be overwhelmed by the background, regardless of the cluster 
size. 

Recall that our estimate of the cluster kernel ax has an un¬ 
certainty of 20% due to the small number of cluster stars. 
Therefore, the concentration of our injected clusters could be 
off by the same amount. We also explore how this uncertainty 
might change our results. The dashed and dotted lines show 
results in cases where our cluster concentration estimate is off 
by 20%. The dashed line shows the result assuming chemi¬ 
cally homogeneous clusters are intrinsically tighter in abun¬ 
dance space by 20%. With a more concentrated signal, the 
signal will have a better contrast over the background. There¬ 
fore, clumps are easier to detect and the probability in Eig. 6 
increases. However, if clusters are more widely spread, they 
are more likely to blend into the background. Therefore, the 
chance of detection decreases, as shown in the dotted line. We 
defer more detail discussions on how this uncertainty changes 
our conclusion to §4.3.2. 

4.3.2. Limits on the CMF 

So far we have only studied the detection probability of an 
individual cluster injected into the APOGEE data. Eor exam¬ 
ple, in the bottom panel of Eig. 6, we derived the probability of 
detecting a cluster as a function of its cluster mass, which we 
will denote as Pdetect(4^ciuster)- In this section, we will propa¬ 
gate this individual cluster statistic to constrain the CME. We 
derive the total number of clusters of different masses from 
Motai and the parameters of the CME. Using this informa¬ 
tion, we can then evaluate the probability of all these pre¬ 
dicted clusters being consistent with the APOGEE observa¬ 
tion, which will then place a limit on the CME. We assume a 
power-law CME with a low-end cutoff of 30 Mq and then con¬ 
strain the power-law slope, a, and the high-end cutoff, Meutoff 
from the comparison with APOGEE data. 

Let’s formulate this idea more rigorously. Given a CME, we 
know that, on average, there are a total ofn = Mtotai/4^ clusters 
spawned where M is the mean cluster mass from the CME. By 
definition, the cluster masses of n clusters follow the CMF, 
which we will denote as Mcluster,i=l,Mcluster,i=2, • • • ,Mcluster,i=n. 


The probability >C(CMF) that all these clusters are consistent 
with the data is 

n 

£(CMF(a, Meutoff)) = n( 1 ^detect (4Tcluster,i)) 5 ( 6 ) 

i=\ 

i.e., none of these clusters exceeds the observation baseline. 
In practice, to save computational time and to ensure a well- 
converged solution, we evaluate the C (CMF) analytically. 

The left panel of Fig. 7 shows the resulting >C(CMF). The 
figure demonstrates that if the CMF slope is shallower than 
a = -2, the APOGEE sample is mostly consistent with a high- 
end cutoff ^ 3 X 10^ Mq (logMeutoff = 7.5). Qualitatively, this 
result should be expected. As shown in the bottom panel of 
Fig. 6, there is a ~ 10% chance that a cluster with ^ 10^ Mq 
will exceed the observation baseline. Recall that if the CMF 
slope a = -2, we have equal contributions from all logarith¬ 
mic mass bins. This implies that the number of clusters with 
mass ^ 10^ Mq is of the order ^ Mtotai/IO^MQ ^ 100. Let 
say there are 50 such clusters, and each cluster only exceeds 
the baseline ^ 10% of the time. The probability that all of 
them would be consistent with the APOGEE observation is 
still extremely unlikely because (90%)^^ <1%. 

This simple illustration also demonstrates two important 
features. First, the detection probability is very low for in¬ 
dividual clusters with masses < 10^ Mq. If the CMF slope is 
steeper than -2, most clusters are not massive. In this case, 
the APOGEE observation provides a very weak constraint on 
the high mass cutoff. As shown in the left panel of Fig. 7, if 
cf < -2, we cannot constrain the CMF cutoff. Even though a 
cluster with mass ^10^ Mq would easily exceed the baseline, 
these clusters are extremely rare if the CMF slope is steeper 
than -2. 

Secondly, since the CMF constraint is derived from the 
product of each detection probability, it is sensitive to ^detect- 
As shown in the bottom panel of Fig. 6, if our estimate of ax 
is off by 20%, it will affect Pdetect^ which in turn could dra¬ 
matically modify our CMF constraint. If clusters are more 
concentrated in abundance space than we have assumed here, 
then that will provide a stronger constraint on each detection 
(dashed line in Fig. 6). Therefore, our constraint on the CMF 
would be conservative. On the other hand, if we have under¬ 
estimated ax by 20%, then the clusters would be more widely 
spread out in abundance space than we have assumed. In this 
case, most clusters would be harder to detect (dotted line in 
Fig. 6). Although not shown, we have checked that, in this 
case, we can only rule out CMF with Meutoff ^ 10^ Mq and 

Nonetheless, independent evidence seems to support our ax 
estimate. Mathematically, the rotation that we performed in 
§3.2 is exactly the same as the principal components analysis 
(see appendix in Ting et al. 2012a). After the rotation, each 
coordinate becomes a principal component of the APOGEE 
abundance space. The APOGEE abundance space has fewer 
independent dimensions than the observed dimensions (An¬ 
drews et al. 2012; Ting et al. 2012a). Some of these 10 
principal components should have very little intrinsic scatter. 
Therefore, some minor axes of the 10-dimensional ellipsoid 
are only due to the measurement scatter. Thus, their spreads 
should be a robust estimate of the measurement uncertainty 
ax. We find that the widths of these minor axes are consistent 
with our ax estimate, showing our estimation of ax is robust. 
Therefore, our conservative CMF constraint is likely to hold. 
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Figure 7. Constraint on the a-enhanced disk CMF. The left panel shows the probability of a CMF being consistent with the APOGEE DR12 data as a function 
of the two CMF parameters, the power-law slope and the upper mass cutoff. The dotted, dashed, and solid lines show the 1-3 sigma limits, respectively. Unless 
the CMF power-law slope is steeper than -2, a cluster mass cutoff >3x10^ Mq is largely ruled out. The right panel shows the maximum cluster mass, Mmax, 
for different CMFs, such that the expected number of clusters Aciuster(> ^max) > 1 (see §4.3.2 for details). We overplot the 1 - 3 sigma limits calculated from the 
left panel. The right panel shows that, in most cases, there is on average less than one cluster with >3x10^ formed in the Milky Way. 


Although not all Mcutoff ^ 3 x 10^ Mq is ruled out, Fig. 7 
shows that a very high Mcutoff is only consistent with the data 
when the CMF slope is steeper than ct = -2. When the CMF 
is steep, the number of massive clusters also decreases pre¬ 
cipitously. Therefore, the cutoff could be very massive, yet 
on average there might be less than one such massive clus¬ 
ter in the Milky Way. A high cutoff does not naturally im¬ 
ply the existence of these clusters. Instead of Mcutoff, Per¬ 
haps a more useful constraint is the maximum cluster mass 
such that we expect to have at least one cluster larger than 
this mass. We denote this maximum mass to be Mmax- As¬ 
suming Mtotai = 6 X IO^Mq, we show M^ax as a function of 
the CMF parameters in the right panel of Fig. 7. As ex¬ 
pected, this panels shows that when the slope is steep, we 
have Muiax ^ Mcutoff, i-^-, the cluster mass cutoff is never 
achieved. When overplotted with the constraints obtained in 
the left panel, the right panel shows that in most cases, only 
clusters with masses < 3 x 10^ could have formed. As the a- 
enhanced disk is believed to form in the first 5 billion years 
(e.g., Haywood et al. 2013), our constraints refer to the por¬ 
tion of the disk that formed at z > 1. 

We have made numerous assumptions in this study, but 
we emphasize that we have always made the conservative 
choices. Therefore, our CMF limit should be robust as long 
as we did not underestimate ax by 20% and the ellipsoidal in¬ 
trinsic distribution is a fair representation of the deconvolved 
distribution. A question remains to be answered: could we 
obtain a significantly stronger constraint on the CMF using 
the current APOGEE data? We would argue that the answer 
is likely no. The bottleneck is intrinsically due to the rela¬ 
tively small number of volume elements in abundance space 
and the low sampling rate. The former is set by the preci¬ 
sion of the abundance measurements and the number of in¬ 
dependent dimensions in abundance space sampled by the 
APOGEE spectra. The APOGEE sampling rate is of the or¬ 
der A^apogee/M totai = 0(10“^), which implies that we would 
only sample one star from a ^ 10^ Mq cluster. The thresh¬ 
old > 10 implies that the minimum cluster mass needed is 
^ 10^ Mq. Therefore, in the most optimistic case, we might 


be able to put a stronger limit by at most an order of magni¬ 
tude. 

How does our stellar CME limit compare to high red- 
shift observations of star forming galaxies? Recent observa¬ 
tions have reported the existence of giant star forming clumps 
within the disks of star forming galaxies at z ^ 2 (Genzel et al. 
2006; Eorster Schreiber et al. 2009; Jones et al. 2010; Liver¬ 
more et al. 2012; Genzel et al. 2013). Gas clumps as high as 
^ 10^ Mq have been observed. There has been much specu¬ 
lation regarding the properties and fate of these giant clumps. 
Some have argued that they are progenitors of globular clus¬ 
ters (Shapiro et al. 2010). They may migrate by dynamical 
friction to the centers of galaxies (Wuyts et al. 2012). It is 
also unclear if the stars forming within these giant clumps 
contain stars that are coeval or share a common metallic- 
ity. What is clear from the results presented in this work is 
that these giant star forming clumps cannot both be mono¬ 
abundance and remain in the a-enhanced disk, at least in the 
portion of the Milky Way observed by APOGEE (i.e., with a 
Galactocentric radius of 3 kpc-15 kpc). Even assuming a total 
star formation efficiency of 1% (simulations and observations 
usually show higher values, e.g. Kennicutt 1998; Elmegreen 
2002; Krumholz & McKee 2005; Evans et al. 2009; Kruijssen 
2012), these gas clumps would have formed clusters that are 
at least 10^ Mq and would have stood out in the chemical tag¬ 
ging experiment presented here if they are chemically homo¬ 
geneous and that they remain in the Milky Way disk, which 
APOGEE is probing. 

4.4. Comparison with previous studies 

The first chemical tagging experiment on dispersed disk 
stars was performed by Mitschang et al. (2014), and the 
tagged groups were subsequently studied in Quillen et al. 
(2015). They studied 714 stars in the solar neighborhood 
from Bensby et al. (2014). Our results agree with their as¬ 
sessments that the identified groups in these studies are prob¬ 
ably not co-natal stars. Each group is unlikely coming from 
a single disrupted cluster, even though the clump members 
might be coeval stars as they share similar abundances. In 
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these studies, the sample includes both ^-enhanced stars and 
low-(T stars, but the sample size is about ten times smaller 
than the APOGEE a-enhanced sample. The sampling rate 
in Bensby et al. (2014) is therefore much smaller than the 
APOGEE (T-enhanced sample. Recall that the sampling rate 
in this study is ^ 0(10“^), and thus, we deduce that the sam¬ 
pling rate in these early studies is ^ 0(10“^)/10 = 0(10“^). 
If groups detected in these studies were to come from indi¬ 
vidually disrupted clusters, the parent cluster would have a 
mass ^ IO^Mq, consistent with the estimates in Quillen et al. 
(2015). 

Simulations from Ting et al. (2015) also disfavor a co-natal 
interpretation of the groups identified in these earlier studies. 
Ting et al. (2015) found that even if such large clusters exist, 
a detected clump in abundance space will still have a sizable 
background component. More importantly, in the case with 
a dominant background, the applicability of previous clump 
search techniques that separate the abundance space into a 
few distinct regions, such as the one proposed in Mitschang 
et al. (2014), or other tree-based methods (e.g. Macfarlane 
et al. 2015) is questionable. Eor those techniques to perform 
well, the background in abundance space has to be negligible 
or first be subtracted. As we have explored in this study, the 
background estimation can be challenging given its complex 
topology in high-dimensional space and the fact that the sig¬ 
nal is usually overwhelmed by the background contaminants. 

5. SUMMARY AND CONCLUSION 

In this study we have exploited the superb APOGEE DR 12 
data, with typical uncertainties of (J^x/fq] and cr[Fe/H] — 0.05 
for 10 elements measured for >10"^ stars, to put in practice a 
first large-scale chemical tagging analysis of the ^-enhanced 
disk. Because the number of stars per 10-dimensional abun¬ 
dance volume is lower in the (T-enhanced disk, we focused on 
that portion of abundance space. 

This analysis required the development of a new, simple 
algorithm for identifying clumps in abundance space, and it 
brought some of the “real life” difficulties of chemical tag¬ 
ging to the fore. Nonetheless, we succeeded in providing 
the first abundance-based constraints on the masses and mass 
functions of chemically homogeneous star clusters in the old 
Galactic disk. 

The methodological steps and results can be summarized as 
follows: 

• We determined and applied a coordinate transformation 
that makes the cluster kernel in abundance space spher¬ 
ical (in 10-dimensions) and have unit variance in each 
dimension. This kernel enables fast error deconvolu¬ 
tions in this transformed abundance space. 

• We generated a model for the intrinsic abundance dis¬ 
tribution of the (T-enhanced disk, presuming it to be 
a highly anisotropic and co-variant ellipsoidal distri¬ 
bution in the above 10-dimensional transformed abun¬ 
dance space. After rotating this coordinate system 
to eliminate the co-variances in this distribution, we 
modeled each dimension independently as the sum of 
two Gaussians. Eitting this to the APOGEE data pro¬ 
vides a first estimate of the shape and volume of the 
error-deconvolved abundance space of the a-enhanced 
Galactic disk. 

• We found that despite the unprecedented quality of the 
APOGEE data, the volume occupied by the stars of the 


(T-enhanced Galactic disk is only ^500 times the vol¬ 
ume of the cluster kernel. Even with abundance un¬ 
certainties of (T[x/fq] and cr[Fe/H] — 0.05 dex, the cluster 
kernel spans > 30% of the abundance width in each ele¬ 
ment abundance dimension. In addition, many of the 10 
element abundances measured by APOGEE and used 
herein are highly co-variant. 

• We developed an algorithm to detect groups of 
chemically homogeneous stars, geared toward the 
background-dominated regime. We found that search¬ 
ing for chemically homogeneous clumps is challenging 
with high backgrounds. The chemical tagging signals 
will most likely come from the peripheral regions in 
abundance space where the background density is rela¬ 
tively low. 

• Using APOGEE data as a detection baseline, we were 
able to constrain the CME in the Galactic ^-enhanced 
disk. We show that this population is unlikely to have 
formed clusters more massive than 3 x 10^ Mq at any 
point in its history. 

Although the current constraints presented in this work are 
limited to very large cluster masses, the results in this paper 
vividly demonstrate the potential of chemical tagging in un¬ 
derstanding the Milky Way properties in the past. With more 
data currently being collected by on-going surveys, we should 
be able to provide much stronger constraints on the CME in 
the near future. 
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